﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace erse253_homework
{
    class Chart
    {

        #region Public Methods

        public void Histogram(double[,] histogramdata, out double[] x, out double[] y, bool cumulative, bool log)
        {

            double min;
            int n;
            double step;

            if (histogramdata != null)
            {
                SetLog(histogramdata, log);
                SetParameters(histogramdata, out x, out min, out n, out step, out y);
                SetXAxis(x, min, n, step);
                SetYAxis(histogramdata, min, step, y);

                if (cumulative)
                    CreateCumulativeDistribution(y);
            
            }
            else 
            {
                x = null;
                y = null;
            }
        }

        public void hScatterplot(double[,] data, double[,] x, double[,] y, int hx, int hy)
        {
            if (data != null)
            {
                for (int i = 0; i < data.GetLength(0) - hx; i++)
                    for (int j = 0; j < data.GetLength(1) - hy; j++)
                        if (i + hx < data.GetLength(0) && j + hy < data.GetLength(1))
                        {
                            if (hx == 0)
                            {
                                x[i, j] = data[i, j];
                                y[i, j] = data[i + hx, j + hy];
                            }
                            else
                            {
                                x[i, j] = data[i, j];
                                y[i, j] = data[i + hx, j + hy];
                            }
                        }
            }
        }


        public void hScatterplot(double[,] data1, double[,] data2, double[,] x, double[,] y, int hx, int hy)
        {
            if (data1 != null && data2 != null)
            {
                for (int i = 0; i < data1.GetLength(0) - hx; i++)
                    for (int j = 0; j < data1.GetLength(1) - hy; j++)
                        if (i + hx < data1.GetLength(0) && j + hy < data1.GetLength(1))
                        {
                            if (hx == 0)
                            {
                                x[i, j] = data1[i, j];
                                y[i, j] = data2[i + hx, j + hy];
                            }
                            else
                            {
                                x[i, j] = data1[i, j];
                                y[i, j] = data2[i + hx, j + hy];
                            }
                        }
            }
        }

        public void Variogram(double[,] data, out double[] x, out double[] y, bool vertical, int tolerance)
        {
            if (data != null)
            {
                SetParametersOutputVector(data, out x, out y, vertical);
                NumericalSolutions solutions = new NumericalSolutions();

                for (int i = 0; i < x.Length; i++)
                {
                    double[,] xh;
                    double[,] yh;
                    GetDataPoints(data, vertical, (i + 1) , out xh, out yh);
                    y[i] = solutions.Moment_of_Inertia(xh, yh, tolerance);
                    x[i] = (i + 1);
                }

            }
            else 
            {
                x = null;
                y = null;
            }

        }

        public void Variogram(double[,] data1, double[,] data2, out double[] x, out double[] y, bool vertical, int tolerance)
        {
            if (data1 != null)
            { 
                SetParametersOutputVector(data1, out x, out y, vertical);
                NumericalSolutions solutions = new NumericalSolutions();

                for (int i = 0; i < x.Length; i++)
                {
                    double[,] xh;
                    double[,] yh;
                    GetDataPoints(data1, data2, vertical, (i + 1), out xh, out yh);
                    y[i] = solutions.Moment_of_Inertia(xh, yh, tolerance);
                    x[i] = (i + 1);
                }
            }
            else
            {
                x = null;
                y = null;
            }

        }

        public void Correleogram(double[,] data, out  double[] x, out double[] y, bool vertical, int tolerance)
        {
            if (data != null)
            {
                SetParametersOutputVector(data, out x, out y, vertical);
                NumericalSolutions solutions = new NumericalSolutions();

                for (int i = 0; i < x.Length; i++)
                {
                    double[,] xh;
                    double[,] yh;

                    GetDataPoints(data, vertical, i + 1, out xh, out yh);

                    y[i] = solutions.CorrelationFunction(xh, yh);
                    x[i] = i + 1;

                }
            }
            else
            {
                x = null;
                y = null;
            }
        }

        public void Covariogram(double[,] data, out  double[] x, out double[] y, bool vertical, int tolerance)
        {
            if (data != null)
            {
                SetParametersOutputVector(data, out x, out y, vertical);
                NumericalSolutions solutions = new NumericalSolutions();

                for (int i = 0; i < x.Length; i++)
                {
                    double[,] xh;
                    double[,] yh;

                    GetDataPoints(data, vertical, i + 1, out xh, out yh);
                    y[i] = solutions.CovarianceFunction(xh, yh);
                    x[i] = i + 1;

                }
            }
            else 
            {
                x = null;
                y = null;
            }
        }

        public double[] Quantile(double[,] data)
        {
            double[] array = SortData(data);

            int n = 200;
            double[] vectorQuantile = new double[n - 1]; ;
            for (int i = 0; i < n - 1; i++)
                vectorQuantile[i] = (array[(array.Length * (i + 1) / n)] + array[(array.Length * (i + 1) / n + 1)]) / 2;

            return vectorQuantile;
        }

        public double[,] CreateSparceMatrix(double[] values, double[] x, double[] y)
        {
            //map to 400 x 500 matrix
            double[,] sparceMatrix = new double[400, 500];
            double m = x.Min();
            double M = y.Min();
            for (int i = 0; i < sparceMatrix.GetLength(0); i++)
                for (int j = 0; j < sparceMatrix.GetLength(1); j++)
                    sparceMatrix[i, j] = double.NaN;

            for (int i = 0; i < values.Length; i++)
                sparceMatrix[(int)(10 * (Math.Round(x[i], 1) - m)), (int)(10 * (Math.Round(y[i],1) - M))] = values[i];

            return sparceMatrix;
        }

        public void CalculateMeanStdDeviation(int size, int divisionsX, int divisionsY, double[,] mean, double[,] stdDevian, double[,] data)
        {
            for (int i = 0; i < divisionsX; i++)
                for (int j = 0; j < divisionsY; j++)
                {
                    CalculateLocalMean(size, mean, data, i, j);
                    CalculateLocalSTDDevian(size, mean, stdDevian, data, i, j);
                }
        }

        public double[,] InflateImage(double[,] data, int size)
        {
            if (data != null)
            {
                int multiplier = data.GetLength(0) > data.GetLength(1) ? size / data.GetLength(0) : size / data.GetLength(1);

                double[,] newData = new double[data.GetLength(0) * multiplier, data.GetLength(1) * multiplier];

                for (int i = 0; i < data.GetLength(0); i++)
                    for (int j = 0; j < data.GetLength(1); j++)
                        for (int k = 0; k < multiplier; k++)
                            for (int l = 0; l < multiplier; l++)
                                newData[i * multiplier + k, j * multiplier + l] = data[i, j];

                return newData;
            }
            else
                return null;
        }

        public void SetParameters(int size, out int divisionsX, out int divisionsY, out double[,] mean, out double[,] stdDevian, double[,] data)
        {
            if (data.GetLength(0) % size == 0)
                divisionsX = (data.GetLength(0) / size);
            else
                divisionsX = (data.GetLength(0) / size + 1);

            if (data.GetLength(1) % size == 0)
                divisionsY = (data.GetLength(1) / size);
            else
                divisionsY = (data.GetLength(1) / size + 1);

            mean = new double[divisionsX, divisionsY];
            stdDevian = new double[divisionsX, divisionsY];

        }
        #endregion

        #region Private Methods

        private static double[] SortData(double[,] data)
        {

            int n = 0;
            for (int i = 0; i < data.GetLength(0); i++)
                for (int j = 0; j < data.GetLength(1); j++)
                    if (!double.IsNaN(data[i, j]))
                        n++;


            double[] array = new double[n];

            for (int i = 0; i < data.GetLength(0); i++)
                for (int j = 0; j < data.GetLength(1); j++)
                    if (!double.IsNaN(data[i, j]))
                    {
                        n--;
                        array[n] = data[i, j];
                      
                    }

            Array.Sort(array);
            return array;
        }

        private static void SetYAxis(double[,] histogramdata, double min, double step, double[] y)
        {
            if (histogramdata != null)
                for (int i = 0; i < histogramdata.GetLength(0); i++)
                    for (int j = 0; j < histogramdata.GetLength(1); j++)
                        y[(int)Math.Round((histogramdata[i, j] - min) / step)]++;
        }

        private static void SetXAxis(double[] x, double min, int n, double step)
        {
            for (int i = 0; i < n; i++)
                x[i] = min + i * step;
        }

        private static void SetParameters(double[,] histogramdata, out double[] x, out double min, out int n, out double step, out double[] y)
        {
            NumericalSolutions solution = new NumericalSolutions();

            if (histogramdata != null)
            {
                min = solution.Min(histogramdata);
                double max = solution.Max(histogramdata);
                n = (int)Math.Round(1 + Math.Log10(histogramdata.Length) / Math.Log10(2));

                //increase the number of bins artificially
                n *= 2;

                step = (max - min) / (n - 1);

                x = new double[n];
                y = new double[n];
            }
            else 
            {
                x = null;
                min = double.NaN;
                n = int.MinValue;
                step = double.NaN;
                y = null;
            }
        }

        private static void CreateCumulativeDistribution(double[] y)
        {
            double[] cy = new double[y.Length];


            for (int i = 1; i < y.Length; i++)
                cy[i] =  y[i];

            y[0] = cy[0];

            for (int i = 1; i < y.Length; i++)
                y[i] = y[i - 1] + cy[i];
        }

        private static void SetLog(double[,] x, bool log)
        {
            if (log)
            {
                NumericalSolutions min = new NumericalSolutions();
                double moveToRight = 0;
                double minx = min.Min(x);
                for (int i = 0; i < x.GetLength(0); i++)
                    for (int j = 0; j < x.GetLength(1); j++)
                        if (x[i, j] <= 0)
                        {
                            moveToRight = -minx + 1;
                            break;
                        }


                for (int i = 0; i < x.GetLength(0); i++)
                    for (int j = 0; j < x.GetLength(1); j++)
                    x[i, j] = Math.Log(x[i, j] + moveToRight);
            }
        }

        private  void GetDataPoints(double[,] data, bool vertical, int i, out double[,] xh, out double[,] yh)
        {

            if (vertical)
            {
                SetParametersOutputMatrix(data, 0, i, out xh, out yh);
                hScatterplot(data, xh, yh, 0, i);
            }
            else
            {
                SetParametersOutputMatrix(data, i, 0, out xh, out yh);
                hScatterplot(data, xh, yh, i, 0);
            }
        }

        private void GetDataPoints(double[,] data1, double[,] data2, bool vertical, int i, out double[,] xh, out double[,] yh)
        {

            if (vertical)
            {
                SetParametersOutputMatrix(data1, 0, i, out xh, out yh);
                hScatterplot(data1, data2, xh, yh, 0, i);
            }
            else
            {
                SetParametersOutputMatrix(data1, i, 0, out xh, out yh);
                hScatterplot(data1, data2, xh, yh, i, 0);
            }
        }
        
        private static void SetParametersOutputMatrix(double[,] data, int hx, int hy, out double[,] x, out double[,] y)
        {

            if (hx == 0 && hy == 0)
            {
                x = new double[data.GetLength(0), data.GetLength(1)];
                y = new double[data.GetLength(0), data.GetLength(1)];
            }
            else if (hx == 0)
            {
                x = new double[data.GetLength(0), data.GetLength(1) - hy];
                y = new double[data.GetLength(0), data.GetLength(1) - hy];
            }
            else
            {
                x = new double[(data.GetLength(0) - hx), data.GetLength(1)];
                y = new double[(data.GetLength(0) - hx), data.GetLength(1)];
            }


        }
        
        private static void SetParametersOutputVector(double[,] data, out double[] x, out  double[] y, bool vertical) 
        {
            if (vertical)
            {
                x = new double[(data.GetLength(1) - 1) ];
                y = new double[(data.GetLength(1) - 1) ];
            }
            else 
            {
                x = new double[(data.GetLength(0) - 1) ];
                y = new double[(data.GetLength(0) - 1) ];
            }
        }

        private static void CalculateLocalSTDDevian(int size, double[,] mean, double[,] stdDevian, double[,] data, int i, int j)
        {
            int numberOfElements = 0;

            if(!double.IsNaN(mean[i, j]))
            {
                for (int k = 0; k < size; k++)
                    for (int l = 0; l < size; l++)
                        if ((i * size + k) < data.GetLength(0) && (j * size + l) < data.GetLength(1) && !Double.IsNaN(data[i * size + k, j * size + l]))
                        {
                            stdDevian[i, j] += Math.Pow((data[i * size  + k, j * size  + l] - mean[i, j]), 2);
                            numberOfElements++;
                        }
            }

            if (numberOfElements > 0)
                stdDevian[i, j] = Math.Sqrt(stdDevian[i, j] / (numberOfElements));
            else
                stdDevian[i, j] = double.NaN;

        }

        private static void CalculateLocalMean(int size, double[,] mean, double[,] data, int i, int j)
        {

            int x = i * size;
            int y = j * size;

            if (i > 0 && i * size + size > mean.GetLength(0))
                x += (i * size + size - mean.GetLength(0));
           
               

            if (j > 0 && j * size + size > mean.GetLength(1))
                y += (j * size + size  - mean.GetLength(1));

            int numberOfElements = 0;
            for (int k = 0; k < size; k++)
                for (int l = 0; l < size; l++)
                    if ((x + k) < data.GetLength(0) && (y  + l) < data.GetLength(1)  && !Double.IsNaN(data[x + k, y + l]))
                    {
                        mean[i, j] += data[x + k, y + l];
                        numberOfElements++;
                    }

            if (numberOfElements > 0)
                mean[i, j] /= (numberOfElements);
            else
                mean[i, j] = double.NaN;

        }

        #endregion
    }
}
